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Abstract 

We demonstrate that the matter density distribution in the surface region is determined well 
by the use of the relatively low-intensity beams that become available at the upcoming radioactive 
beam facilities. Following the method used in the analyses of electron scattering, we examine how 
well the density distribution is determined in a model-independent way by generating pseudo 
data and by carefully applying statistical and systematic error analyses. We also study how the 
determination becomes deteriorated in the central region of the density, as the quality of data 
decreases. Determination of the density distributions of neutron-rich nuclei is performed by fixing 
parameters in the basis functions to the neighboring stable nuclei. The procedure allows that the 
knowledge of the density distributions of stable nuclei assists to strengthen the determination of 
their unstable isotopes. 
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1 Introduction 


Experimental facilities for unstable nuclei are receiving much attention in the recent years. Radioac¬ 
tive Ion Beam Factory (RIBF) is under construction at RIKEN in Japan, and Rare Isotope Accelerator 
(RIA) facility is being proposed in USA. The facilities are expected to provide much data which are 
expected to reveal interesting properties and dynamics of unstable nuclei 1^, ^. Among various prop¬ 
erties, one property of interest is the one-body matter density distribution. The density distribution is 
a fundamental quantity of nuclei and serves as an important measure to test how well we understand 
nuclear structure ||] , and will be experimentally investigated at an early stage after completion of the 
facilities. It is perhaps an appropriate time to make a close investigation how well one could determine 
the one-body matter density distributions from the data expected to emerge from the facilities. As a 
concrete example, we examine the type of the experiments under consideration at the RIBF. 

Unstable nuclei of neutron-rich side (neutron-rich nuclei) are often characterized as those with a 
large surface region generated by loosely bound valance nucleons iiiH- We will give a special 
attention to the determination of the matter density distributions in the surface region, though we 
will have to treat, of course, the entire region of nuclei. The surface is the region where p(r) /p(0) 
varies from 0.9 to 0.1 if we have a density distribution of monotonically decreasing function in r, such 
as the Woods-Saxon density, in mind Q. More elaborated definitions of the surface can be found in 
the literature Q. 

In this work, we focus on the question, how the precision in the determination of the matter den¬ 
sity distributions of neutron-rich nuclei is related to the accuracy of experimental data, specifically 
of proton-nucleus elastic scattering differential cross sections. For the purposes of making the exam¬ 
ination definitive, we perform a model calculation, by artificially generating data from given density 
distributions and then by making their statistical analysis. 

The question which we address is an old problem, and the basic method that we use has been well 
developed for extracting nuclear charge density distributions from high-energy electron scatterings 
and muonic atom data over years |l^, |^, |^, We make no new contribution to the method, 
but apply it to the determination of the nuclear matter density distributions, focusing on the surface 
region of neutron-rich nuclei. As we discuss in the next section, the method has been applied to the 
determination of the nuclear matter density distributions but with less rigor, and no use has been 
made to the case of neutron-rich nuclei. The present work is novel in the careful examination of how 
the method is applied to the nuclear-matter-density-distribution problem. 

In this work, we follow the traditional approach of the least-square method |Q and do not discuss 
other approaches often referred to as Bayesian methods, such as maximum entropy technique. These 
approaches involve issues of subjectivity |^, |^, and require a separate study. We are currently 
investigating their relevance to our problem. 

The contents of this paper are as follows: We describe our approach in comparison to previous works 
in Section |^, and the formalism used in our analysis in Section There, we introduce pseudo data 

Numerical results are presented in Section J4 


3.2. 


in Subsection 3.2 and error estimate in Subsection 

including discussions of pseudo data in Subsection 1^ and those of error estimate in Subsections I 
and |4.3| . We summarize our results in Section ||. Appendices include a brief discussion of basis functions 
in Appendix ^ a list of the expressions for the density and the data in Appendix ^ and a summary 
of formula and numerical results concerning the Fourier-Bessel basis functions in Appendix 

Some preliminary results of this work have been previously reported p^ . This paper provides a 
full account of our work 0. The major added contents are the discussion of the completeness error, 
the way how to generate pseudo data, and the other confirmation of the fitted results with the use of 
the x^-distribution. We have also added several Appendices to the previous report. 


^This paper is the revised version of the preprint, RIKEN-AF-NP-413, Nov., 2001, with the same authors. 
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2 Overview of Our Approach 


In this investigation, the density distributions are determined in a “model-independent” way. As the 
model-independent determination is the central feature of our procedure, we wish to clarify at the 
outset what we mean by the model-independent procedure. It is model-independent, because the 
density distribution is determined as a linear combination of basis functions. We do this by avoiding 
the use of specific analytic forms of the density distribution, which is the common practice in this 
field 0, 0 [T|, [T| 


We will explicitly assume, however, that the density distribution is directly related to the differ¬ 
ential cross sections through the eikonal approximation in the lowest-order optical potential, so-called 
t-p form fl^ . We assume this because of convenience and simplicity. The assumption is basically 
independent of our semi-quantitative conclusion on the question of how the density determination 
depends on experimental accuracy. That is, we expect that the same conclusion would emerge if we 
were to adopt a more elaborate relation. We emphasize that the issue we wish to address in this paper 
is this question of the dependence of the density determination on the experimental aceuracy, not on 
the relation between the density and the cross section themselves. 

The elastic scattering cross section data that we will examine in our model calculation are gener¬ 
ated artificially, and will be referred to as “pseudo data” [1^ 1^ At first, they are generated 


from a specific density distribution (of the three-parameter Fermi form, Eq. (4.1)) by the use of the 
above described procedure with the eikonal approximation and are then shuffled statistically at each 
data point with an assigned experimental uncertainty (corresponds to a Gaussian width), simulating 
the situations expected in future experiments at the RIBF. We will prepare 25 sets of such shuffled 
data, each of which yields the best-fit density distribution with its uncertainty. The originally assumed 
density distribution is taken to be the true density distribution. The average of the density distribu¬ 
tions determined will be close to the true distribution, and their deviations from the average will give 
us the statistical uncertainty associated with the model-independent determination |Q. The differ¬ 
ence between the true and averaged distributions corresponds to the systematic uncertainty Q. The 
systematic uncertainty corresponds to the limitation of the model-independent procedure associated 
with our (model-independent) representation of the density distribution. The systematic uncertainty 
could also include the uncertainty associated with the relation between the density distribution and 
the differential cross sections. In this work, we will not examine this part of the systematic uncertainty 
by applying the same relation for generating the data and fitting to them. 

We apply a particular combination of Gaussian functions |pi| as a model-independent represen¬ 
tation of the density distribution. We do this for simplicity and for stability of computation, but 
a particular choice of the basis should not alter the final conclusion. That is, the nuclear density 
distribution determined at the end should fall in within the systematic error that is noted above. We 
will confirm this by comparing the results of our Gaussian basis functions with those of the more 
established Fourier-Bessel basis [ p^ . 

There are several well-known basis functions. The expansion in terms of ^-functions was proposed 
originally by Lenz [^, and extended by Lenz and Friedrich |^. Since the early 70’s, Sick and 
his collaborators have been analyzing the electron scattering data using another set of Gaussian 
basis functions, the expansion in a sum of Gaussians (SOG) Q, with great success ||^, p^ . 

These Gaussian basis functions are a series of Gaussians of multi-center. With their analyses, they 
demonstrated that the determinations of the charge density distributions by the phase-shift analyses 
is reliable and accurate by comparing the electron-nucleus with the positron-nucleus elastic scattering 
[psi , and the charge density distributions of ^°®Pb was determined with an uncertainty of the order of 
Ap{r)/p{r) = ±1 [%] Ip^. Recenty, Burov et al. analyzed the electron scattering data by expanding 
the density distribution in terms of the symmetric Fermi distribution |2^ and its derivative series . 

The three-parameter Fermi distribution that we use mimics the matter density distribution of ^®Ni 
and the relativistic mean field result mimics that of ^®Ni, which we chose as examples of stable and 
unstable nuclei, respectively. We assume here that the matter density distribution of ^®Ni behaves 
similarly to its charge density. From the intense studies of the charge density distribution |^, po| , 
it was found that the charge density of ®®Ni shows the less oscillating structure than the predictions 
of Hartree-Fock calculations pffl . This validates our choice of this simple functional form (See also 
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in Appendix The pseudo data are generated for the proton incident energy of 1.047 [GeV]. The 
energy is chosen to be the same as that of the proton-nucleus scattering experiments at SATURN, 
Saclay ^ because the experimental data and their analyses have shown that the eikonal 
approximation yields reasonably realistic results (See Fig. 1 in the next section). 

In this work, we assume that the proton-nucleus elastic scattering differential cross sections and 
the nuclear density distribution are exactly related by the eikonal approximation with the t-p form. 
With this simple form, one can see analytically that the nuclear radius determines the oscillation 
and the diffuseness does the exponential decrease with momentum transfer of the elastic cross section 
[Q. This helps one to obtain the insight to the relation between the quality of data and the fitted 
density distributions. Of course, the eikonal approximation is an approximation, and the degree of its 
validity is sometimes questioned, especially when it is used with the t-p form. It can be numerically 
improved by the use of the partial-wave decomposition method, and the t-p form can be also improved 
by incorporating higher-order effects into the potential. Our objective of finding the relation between 
the precision of the density determination and the experimental accuracy certainly depends on how 
reliably the density and the data are related each other. For the sake of making our examination as 
clear as possible, however, we isolate the issue of the precision-accuracy in this model calculation, by 
artificially assuming the density-data relation to be exact. In realistic model-independent analyses of 
actual experimental data, we must address the latter issue and take account of the reliability of the 
relation of the two quantities. 

Once the density-data relation is assumed to be exact, the model-independent analysis can provide 
a definitive result for the precision-accuracy issue at any energy. The result will become, however, 
unrealistic as we go down in the energy. For a model calculation at lower energies, such as those 
around 500 [MeV] where various experiments are being planned at the RIBF, it would be desirable to 
use more rigorous relation based on the partial-wave decomposition and an improved optical potential. 
We are planning to extend the present work to these energies by the use of partial-wave decomposition. 

There have been efforts on the determination of the nuclear (matter-) density distributions 0 
0 0 ^1, which are, of course, founded on the charge density studies. As reviewed by Alkhazov 
et al. [^1 and by Chaumeaux et al. in 1978, and by Batty et al. 0 in 1989, some works 
include model-independent analyses of the experimental data: after the first attempt of the model- 
independent approach to the matter density distributions by Brissaud et al. [Q , Ray et al. ^ , 
Hoffmann et al. [^, and Alkhazov et al. |^2| studied the matter density distributions along similar 
lines. More recent work is by Starodubsky et al. Gils et al. analyzed the data of the a-nucleus 
elastic scattering, and obtained the quite similar surface density distribution as those from the proton- 
nucleus scattering data within the statistical uncertainties 0]. None of these works has addressed, 
however, the issue of determination of the surface density distribution associated with nuclei near 
the drip lines. The recent analysis by Alkhazov et al. |Q of the GSI experiment of helium isotope 
addresses this issue, but it is a model-dependent analysis and its conclusion is in question under the 
light of a model-independent analysis . 

Our work is the first model calculation that examines closely the information content of the ex¬ 
perimental data, including the question of the systematic error as well as the statistical error, in the 
model-independent approach. 

Although they are model-dependent, the matter density distributions of neutron-rich nuclei have 
been deduced recently, e.g. from the low-energy proton elastic and inelastic scattering 46, and 
from the interaction cross sections Eq, 


50 
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3 Formalism 


3.1 Fitting Procedure 

For the sake of notations and completeness, we outline our formalism in this section. It is basically 
the same as that in the review article by Friar and Negele ||^. The formalism itself was studied also 
in detail by Dreher et al. [ p^ . 

We approximately express the one-body matter density distribution that is assumed to be true, 
/Otrue(i’), in M-terms of basis functions, /n(i‘), (n = 1, • • •, M): 

M 

Ptrue(l‘) — Pm{^) = Cn fni^), (3-1) 

n—1 

where (r) is the fitted density distribution. Here, as we do not use a particular analytical functional 
form of the density, our approach is model-independent. Strictly speaking, there is no completely 
model-independent approach, because the basis functions include parameters to be chosen and only 
a finite number (M) of the basis functions can be used. The systematic errors are generated due to 
this incomplete model-independence of the approach. We will describe this point fully later. 

For simplicity, we assume that the nuclei are not deformed and examine only the radial dependence 
of their distributions. The deformation can be treated in a similar model-independent way, but will 
make the procedure more complicated. The fitted density is normalized to the mass number. A, as 

J dr PMir) = 47r^ r^dr p^{r) = A. (3.2) 

The expansion coefficients, {Cn}, are determined by minimizing 


X^aCn}) 


OL—1 ^ ^ 


da] 


fit 




■{{Cu}) 


(3.3) 


under the variations of the finite set {C„}. Here, dua/dH’s are the differential cross sections of proton- 
nucleus elastic scattering at the center-of-mass angle, Oa- IV^ata is the number of data points, and Cq, 
is the uncertainty or error associated with the data at 9^- da'^^*'^/dil's are the pseudo data that we 
generate, and da^^/dfl{{Cn}ys are those calculated from the fitted density distribution. We carry 
out all least-square minimization using MINUIT |5^ in the CERN Library. 

The dependence of da^^/dfl on {Cn} involves complicated physics. In this work, we use the 
first-order optical potential in the Glauber approximation, which establishes the relation between the 
nuclear density and the differential cross section 0- The relation is well known, but we list some of 
basic formulae in the relation in Appendix ^for completeness and clarification. As noted in Section 
we assume the relation to be exact and generate pseudo data using it. No real experimental differential 
cross section datum enters our calculation and its analysis. As such, it is logically immaterial how 
reliable the relation itself is, but it becomes an issue how realistic or relevant to future experiments. 

Figure I compares the real experimental data of proton-®®Ni elastic scattering cross sections at 
1.047 [GeV] from SATURN, Saclay ||^, Q with the calculation. Note that in our later model 
calculations, we will use the same nuclear density at the same kinetic energy as used here. Agreement 
between the data and the calculation is not too great, but is not too much off either. The agreement 
could be improved by including the Coulomb interaction and also various higher-order corrections, 
such as the center-of-mass correction and nuclear correlation 1^, |^. The improvement would reduce 
the unpleasant sharp dips in the differential cross sections [pS and is certainly desirable, but again, 
is not essential to our model calculations. 

As the basis functions of Eq. (^^), we use Gaussian basis functions with a cut off at R: 

fn{r) =exp{-r^/rl) e{R-r), (n = l,---,M), (3.4) 

where r„ is the size parameter, defined in the geometrical progression, rn = (n = I, • • •, M) 

pl| . These Gaussian basis functions differ from those used previously by Sick , and will be referred 
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to as “Kamimura-Gauss basis functions” after the author who used them extensively H- We lists 

the basic formulae of the Kamimura-Gauss basis functions in Appendix ^ _ 

In Subsection 4.1 we will show that the completeness error defined in Subsection 3.3 of the 
Kamimura-Gauss basis functions is very small, which is the same order of magnitude as that of 
the Fourier-Bessel ones shown in Appendix |^. 


3.2 Pseudo data 

We create the total 25 sets of the pseudo data. We choose the number of 25 as the minimum number 
for statistical discussions. Each set mimics a series of experimental data obtained by an experiment, 
and consists of the proton-nucleus elastic scattering cross sections at the center-of-mass scattering 
angles that are common to all pseudo-data set. No Coulomb interaction is included, for simplicity. 

The application of the pseudo data to estimate the density distributions is not new 1^, In 
this work we study how far we can access the matter density distributions with varying the nuclear 
beam intensity, having in mind the neutron-rich nuclei. This is new. 

At each scattering angle, 6a, we thus generate 25 cross sections. They are generated artificially, 
but as realistically as possible. The following steps are taken to achieve this objective. First, we 
calculate the “true” cross section, da^'^^/dV,, from the given nuclear density of a three-parameter 
Fermi distribution form, as described in Appendix Second, a set of 25 cross sections are generated, 
randomly distributed about each “true” cross section. Third, an uncertainty is assigned to each datum 
generated. 

Let us elaborate on these steps. In the second step, 25 cross sections are generated following the 
Gaussian distribution, 


r/ 2 ^ 1 

f{x-n,a ) = -—- exp{ —^-} 




2a2 


where 


d- = 


dal 


dn ’ 

d’i \J Ayield(^Q 


Here, we define the yield count, 


Ayield(^a) = 


No da^r 

Bc] (ill 


(3.5) 

(3.6) 

(3.7) 

(3.8) 


where Bq/Nq [mb/str] is the yield parameter, which is the magnitude of the cross section that yields 
a unit count in the measurement. We set Nq = 10, so that 


a = 0.316^, when 


dcrf, 




= Bo 


(3.9) 


and we examine three cases of Bq = 0.1, 1.0, and 10.0 [mb/str]. Note that the datum less than Bq/Nq 
is rejected, because A^yieid(0a) of Eq. (3.8) becomes less than one. 

I n th e third step, we assign an uncertainty to each pseudo datum, using the expression of a in 
Eq. (^) with Eq. (^), except for da//™/dVl being replaced by the pseudo datum itself. 

The center-of-mass scattering angle is taken to be greater than 4°, based on the expected exper¬ 
imental setups (of the inverse kinematics) at the RIBF. The set of the angles, {6a}, is chosen with 
1.0° separations, but the separation is reduced near the forward direction as described in Appendix ^ 
The number of the angles (and thus the number of the data in each pseudo-data set, iVdata) is then 
made to be more than 20, which provides the number of degrees of freedom, Ndof (= A^data — Af), 
suitable for the least-square fitting procedure that we will carry out for M = 10. 


The pseudo data are generated by following the above recipe, as will be described in Subsection 4.2 
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3.3 Errors and Uncertainties 

In this subsection we discuss how to estimate the errors involved in the determination of density 
distributions. There are three kinds of errors; the completeness error, the systematic error, and the 
statistical error |]^ . 


3.3.1 Completeness Error 

The completeness error is a measure of how closely the density distribution is reconstructed by the 
use of a finite number of the basis functions |l^. We introduce the relative completeness error as 

= |ptrue(?’) - PM{r)\ /ptrue(?') , (3.10) 


where 


M 

PM{r) = '^Cn fn{r). 

n—1 


(3.11) 


{Cn} is determined to minimize (ptrue(r) — pM{r)Y under the variations of Cn’s. Note that fn(rYs 
are linearly independent set, but are not necessarily orthogonal among each other. Since we use a 
finite number of the basis functions, the completeness of the set is not the issue here. 

{(7„} is given by 


as the solution of 


where 


a 


M 

n—1 


M 


Fn^mCm — Un, 




Vn 



r’'dr fn{r)fm.{r), 



r’'dr fm{r) ptme{r). 


(3.12) 


(3.13) 


(3.14) 

(3.15) 


Here, I is the phase-volume parameter to specify the relationship among the basis functions, as de¬ 
scribed in Appendix In this work, we choose I = 0. 

The relative completeness error introduced here will provide a useful means for assessing how well 
the surface distribution is determined. In this work, we use the Kamimura-Gauss basis functions for 

Unir)}. 


3.3.2 Statistical Uncertainties and Systematic Errors 

When the true value, Wrue, is known, the mean square error of a measured value, Aldataj is defined 
as 0, 

((Xdata - Wtrue)') = ((^data " (Wdata))") + ((^data) " ^true)", ( 3 . 16 ) 

= cr 2 -k 6 ^ ( 3 . 17 ) 


where the bracket, (A), implies the sample mean of A. a is the statistical error, and h is the systematic 
error or the bias. Generally, the error contained in data is an independent sum of the statistical and 
systematic errors. The statistical error shows how large the fitted results fluctuate around the sample 
mean, while the systematic error shows the deviation of the sample means from the true value J 
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In our case, the statistical error of the fitted cross sections at a 0cm is 


and the systematic error is 




b\{Pcn\) = 


_ / da 


fit 


,da 


fit 


' \ dfl 




,da 


fit 


dcr‘ 


dn dn 

The sample mean denoted as (• • •) in Eqs. (3.IS) and ( 3.19| ) at 0cm is explicitly 

.da^\ = 1 A 

n 




n ^—' dil ’ 

i—1 


where the number of the pseudo-data set n is 25. 

We define the statistical and systematic errors for the fitted density distributions in the 
way. 


^p(rf = ((pmC?-) - (PmW))\ 

bp{rf = ((PmW)- 0true(?'))^ 
where the sample mean is defined in the same way as in Eq. (|3.20D. 


(3.18) 


(3.19) 


(3.20) 

similar 

(3.21) 

(3.22) 
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4 Numerical Calculations and Discussions 


4.1 Completeness Error 


We first examine the completeness error, Eq. (3.10), of the Kamimura-Gauss basis functions, taking 
®®Ni density as an example. The true density distribution, ptrue(r), is taken to be of the three- 
parameter Fermi distribution form |p4| , 


Pir) = po 


1 -I- w{r/c)'^ 


1 -I- exp{(r — c)/a}’ 


(4.1) 


where w = —0.1308, c = 4.3092 [fm], and a = 0.5169 [fm]. po is the normalization constant determined 
by 


J dr p{r) = Itt y r'^dr p{r) = A. 


(4.2) 


Using Eq. (El above as ptrue(r) in Eqs . (|3.12D - ( |3.15D , we can determine {CVi}, or PM(r) from 
Eq. (3.11), and thus A/9^™(r) using Eq. ( 3.101 ). Our basis functions, the Kamimura-Gauss basis 
functions, however, depend on the size parameters, ri and tm, defined in Eq. (3.4). All quantities 
involved for determining A/9^™(r), such as {Cn}, thus also depend on these parameters. We determine 
the values of ri and tm for a fixed M, so as to also minimize 

{ptrue{r) - pM(r)}^, 


(4.3) 


for the variations of ri and tm- This minimization is carried out numerically. 

Table 1 lists the best-fit values of ri and tm for different M, and Fig. 2a) shows the corresponding 
relative completeness error, Ap^“(r). The radial range for the fit is i? = 10 [fm]. In this completeness 
calculation, we impose no normalization condition on pM^r), and thus the integrated value of pm^t) 
also serves as a measure of goodness of pM{r)- The integrated value. Am, is listed also in Table 1. 
Note that the exact integrated value is 58 in this case. 

Judging from the Am value in Table 1 and Ap\°’^{r) in Fig. 2a), we see that the appropriate M 
would be greater than 10. As we will see in Subsecion [4.2| , the uncertainty argument gives M to be 
about 10. We will choose M = 10 in the rest of this work. 

Figure 2a) shows that Ap'^^ir) increases in average as r increases. For M = 10, Ap^™(r) increases 
to about 0.01 around 7 [fm] from 0.0001 at the central region. Though the increase is appreciable, it 
would not be alarming as the density is down by the factor of about 1000 around 7 [fm], compared 
to the central region of the density. This behavior is quite similar to that of the Fourier-Bessel basis 
functions (see Appendix 

We observe that while the ri value depends strongly on M, the tm rather weakly. Because of this 
strong dependence, one may wonder how precisely the size parameters have to be chosen while keeping 
the completeness error to be reasonably small. In order to address this issue, we choose a fixed set 
of the size parameters, ri = 0.5 [fm] and tm = 6.0 [fm], which deviate from the best-fit values for all 
M shown in Table 1. Figure 2 b) illustrates Ap^™(r) for this choice of the size parameters. Though 
Ap^™(r)’s in this case are larger than those for the best-fit parameter values, we consider Ap^™(r) 
in the order of 0.01 to be tolerable. Am comes out to be quite close to the exact value, off by 3 [%] 
even in the worst case of M = 5. 

The mild dependence of the completeness error on the size parameters will be further demon¬ 
strated as we a pply the same best-fit size parameters for ^®Ni to the case of ^®Ni. As described in 
Subsection iA, we can exploit this mild dependence when we extend the model-independent analysis 
from the relatively well-known stable nuclei to poorly-known neutron-rich nuclei by keeping the same 
size parameters. 

The Kamimura-Gauss basis functions are not orthogonal among each other. Generally, however, 
the minimum point of the least-square does not change in the M dimensional space of the basis 
functions under a rotation of the space coordinates attached to the basis functions. The minimum 
point thus remains the same after the Kamimura-Gauss basis functions are orthogonalized by a linear 
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M 

ri [fm] 

Tm [fm] 

Am 

5 

0.6 

5.8 

56.255 

10 

1.0 

5.5 

58.030 

15 

0.7 

5.6 

57.999 

20 

0.3 

5.7 

58.000 


Table 1: Dependence of the best-fit parameters of ri and tm in the Kamimura-Gauss basis functions 
and of the integrated value of the density Am on the number of the basis functions M. The radial 
range of the fit is i? = 10 [fm]. 


Group 

Bq [mb/str] 

Cmin [%] 

A0„iin [mrad] 

Rate [sec 

A 

0.1 

2.0 

2 X 0.3 


B 

1.0 

3.3 

2 X 0.3 

10^ 

C 

10.0 

6.7 

2 X 0.3 

10 


Table 2: Parameters of the pseudo-data-set group A, B, and C, respectively. Bq is the magnitude of 
cross section when the yield count is Nq. Nq = 10 in this work. Cmin [%] is the minimum uncertainty. 
Akinin [mrad] is the minimum step. 


transformation, say, through Schmidt orthogonalization method. Note that this is the case even if a 
nonlinear transformation were to be used. The error matrix at the minimum point does change, of 
course, under a transformation, linear or otherwise, by altering the correlated nature of the uncertain¬ 
ties. In the preceding discussions, we have addressed to this issue by performing numerical calculations 
by artificially moving away from the minimum point. This approach would be more reasonable, as 
we do not know the true distribution and thus the exact location of the minimum point in practical 
applications, neither. 

Before closing this subsection, we show the dependence of the cross sections on M in Fig. 3. The 
Kamimura-Gauss basis functions give satisfactory convergence for M > 10. To go beyond 20° it gives 
nice diffraction patterns, especially for M > 15. The solid curve (M = 20) is on another dotted curve 
which is drawn by using the original ptrue{r) for 0cm < 20°, and slightly differs for 0cm > 20°. This 
behavior is better than that of the Fourier-Bessel basis functions, (see Appendix ^). 


4.2 Pseudo Data 


As we described in Subsection 3.2, we generate 25 pseudo data at each scattering angle. The values 


of the pseudo data change when the yield parameter Bq is changed. We carry out our analysis for 
three different values of Bq = 0.1, 1.0, and 10 [mb/str]. Each value of Bq [mb/str] provides different 
sets of the pseudo data, each set consisting of 25 data at each scattering angle. We denote the set for 
Bq — 0.1 [mb/str] to be the pseudo-data-set group A, the set for Bq = 1.0 [mb/str] to be the group 
B, and the set for Bq = 10.0 [mb/str] to be the group G. We take the minimum uncertainty Cmin = 
2.0, 3.3, and 6.7 [%] for pseudo-data-set group A, B and G, respectively (Table 2). The definitions of 
these quantities are given in Appendix |^. 

As an example, we plot a pseudo-data set of the pseudo-data-set group A in Fig. 4. The solid 
curve is the same as in Fig. 1. The kinematical conditions are: 


0max<25° O gmax < 750[MeV]. (4.4) 

A0~1.O° Ag~20[MeV], (4.5) 

where Ag is an increment in momentum transfer. Since we take a step in 0 smaller than 1° in the 
forward direction, Ag can be smaller. 

Some of the parameters of the expansion of the density distribution are related to the kinematical 
conditions of the experiment [|^ |I^, H^, . The followings are the rough estimates for them: rmax 

is the maximum distance to probe by the experiment, and r^ax — 7r/Ag ~ 30 [fm]. We take R 
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= 10 [fm] (< Tmax)- Since we are interested in the surface region, not in the tail, this value of R 
is enough to cover the region. This choice is consistent with the charge density distribution studies 
for the medium-heavy nuclei, 50.52.54 qj. resolution of the fitted density distribution is Ar 

= 7r/(7„iax — 0.85 [fm]. The dimension of the model space is M ~ qmaxR/T^- We estimate M = 10 ~ 15. 

In order to make our discussions more realistic, we compare different pseudo-data-set groups with 
experimental setups expected at the RIBF. We first introduce two quantities, the counting rate, J(0 q) 
[sec“^], and the time for measurements, T [sec]. The multiplication of the two gives the yield count, 

NyieidiOc.) = j(0„) T (4.6) 

in Eq. (^). J(0 q) is expressed as 

J true 

J{0o) = War Jbearrr (4.7) 

where War is the number of the particles in the reaction region of the target, jbeam [cm“^ sec“^] is 
the current of incident particles, /dH. [cm^ str“^] is the cross section at 9a, and Afl [str] is the 

solid angle covered by the detector. 

At the RIBF, experiments of proton-nucleus elastic scattering are performed in the inverse kine¬ 
matics that the incident nuclear beam collides with protons in the target in the laboratory frame. 
Consider a CH 2 target with the thickness of 10 [mg cm“^]. In terms of the Avogadro number Na 
= 6.03 X 10^^ [mol“^], the number of the protons in the target is War = 2 x Na x(10 x 10“^/I4) 
= 2 X 0.43 X 10^^ [cm“^]. Let us assume the radius of the incident beam to be 1 [cm] and the number 
of the detectors to be 10 covering 50 [mstr] each (that is, AH = 0.5 [str]), and also take /did 

= 0.1 [mb/str], and the rate of particle production to be 10^ [sec“^] (that is, jbeam = IO^/tt [cm“^ 
sec“^]). J{6a) is then in this case 

103 

J{9a) = tt [cm^] X 2 X 0.43 x 10^^ [cm“^] x - [cm“^sec“^] 

TT 

X 10“^® [cm^str“^] x 0.5 [str] 

= 4.3 X 10"® [sec“^]. (4.8) 

If we suppose 1-week measurement, i.e., T = 6.048 x 10® [sec], we have Wieid(^a) = 26.0, or in the 
order of 10 counts. Therefore, 0.1 [mb/str] roughly corresponds to a 10 yield for the above conditions, 
which is pseudo-data-set group A. The correspondence of the rate of particle production to various 
pseudo-data-set groups is summarized in the Table 2. 

Note that the expected production rate of ™Ni at the RIBF is about 10® [sec“^], while that of 
™Ni is about 10 [sec“®]. The former corresponds to the pseudo-data-set group A and the latter to the 
group C. These correspondences are still valid if the target thickness becomes 5 [mg cm“^], because 
we have estimated them by considering Ayieid(6*a) = 26.0 to be the order of 10 counts. 


4.3 Least-Square fitting 


For the three pseudo-data-set groups of proton-®®Ni elastic scattering generated in the previous subsec¬ 
tion, we perform least-square fittings of Eq. (|3.3D in Subsection 3.1 to obtain the density distributions. 


p{r). In the fittings, the normalization of p(r) is treated as one of the pseudo data with an uncertainty 
of about 0.2 [%], far smaller than the uncertainties assigned to the other pseudo data. Note that the 
normalization of the basis functions is automatically taken care of in this way, because such a high 
precision datum works as a constraint. The numerical minimization is done by the use of MINUIT 


There is another usage of the normalization. Without putting a high precision number to the 
normalization, one may make the fitting first, and calculate the normalization after the fitting is 
completed. This normalization should fluctuate around the true number, and the degree of this 
fluctuation could work as a measure of the goodness of the fitting. We do not take this way, because 
we are interested in the statistical and systematic errors of the density distribution for the fixed 
normalization. 
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The numerical results for the pseudo-data-set group A, B, and C are shown in Figs. 5, 6, and 7, 
respectively, for i? = 10 [fm] and M = 10. As the particle production rate decreases from Fig. 5 to 
Fig. 7 (see Table 2), quality of the pseudo data becomes worse accordingly. Each hgure consists of 

a) 25 htted cross sections (shown with solid curves) with all pseudo data (crosses with bars), and 

b) 25 htted density distributions, (solid curves), with the original density distribution (the dotted 
curve). The density distributions are shown for Pm(^) including the radial phase space of . It 
is this form that the data actually constrain. In the hgures we show the statistical and systematic 
errors divided by the true density distributions. One can see that the relative uncertainties of the 
htted density around the true density itself without from those hgures. The average values of the 
least square divided by the number of degrees of freedom, y^/DOF, and of the RMS radius, (r^)^/^, 
obtained from the 25 diherent pseudo-data sets are also shown in the hgures. The number of DOF is 
(= iVdata -M). 

For all pseudo-data-set groups, the density distributions, have larger uncertainties in 

the central region than in the surface region. This result is in agreement with the previous works 
@1, © and is analogous to what was found in the determination of charge density 

distributions (for example, ). The original distribution drawn by the dotted curve is included in 

the error-band, which signihes that our procedure consistently provide probable density distributions 
within the uncertainty. Even with a pseudo-data-set group of lower yields, it seems that the surface has 
smaller uncertainty than the center. It is reasonable that the surface region of P%{t) is determined, 
because the cross section carries the information of this region 

In each case, y^/DOF is about unity. For the pseudo-data-set group B we have checked that the 
least square, y^, of each pseudo data set distributes according to the y^-distribution (Appendix 
This confirms that our numerical simulations make sense statistically. The of each pseudo- 

data-set group has less than 2 % uncertainty, including the original value of 3.764 [fm]. 

Comparing Figs. 5, 6, and 7, we find that while the differential cross sections fluctuate for 0cm > 
20°, the nuclear densities fluctuate for r < 2[fm], and that the fluctuations become more appreciable 
again as the quality of the pseudo data decreases. Note that when the number of terms of the 
Kamimura-Gauss basis functions, M, is increased, the y^/DOF becomes closer to unity but the 
fluctuations increase, as similar observation was made previously |]l3| , . 

In order to see explicitly the role of large 0cm or large-g pseudo data, we remove the pseudo data 
of 0cm > 15° from the pseudo-data-set group C, and call the truncated data set to be the group C’. 
Figures 8 a) and b) show that the fluctuations by the pseudo-data-set group C’ are greater than those 
by the group C both in the cross sections and in the density. Note that the mean value and the 
variance of changes much less, because the density in the surface region mainly affects (r^)^/^. 

The statistical uncertainty and the systematic error are examined in Figs. 9 and 10. We see 
that as the quality of the pseudo data decreases, the statistical uncertainty tends to increase but the 
systematic error does less. Furthermore, Fig. 9 shows that both of the uncertainty and the error in 
the surface region of our interest, r = 3 ~ 5 [fm], are smaller than in the other regions, which is in 
consistent with the previous works ]^ , This is one of our major results. 

The quantities examined here are the relative quantities. One should not be alarmed with the 
somewhat large values for r > 6 [fm], because they are associated with quite small values of the 
density. A similar caution is needed with the cusps appearing in Fig. 10, because they are associated 
with the dips in the cross sections. The cusps are expected to reduced once the Coulomb interaction 
is included in the analysis. 

4.4 Neutron-rich nuclei 

We now examine an extension of our analysis to the case of neutron-rich nuclei by taking the case 
of the proton-’^®Ni elastic scattering as an example. Here, we use the pseudo-data-set group C to 
generate its pseudo data, as the expected particle production rate of "^^Ni at the RIBF, 10 [sec“^], 
corresponds to the case of the group C. As the true density, we take the vector part of the relativistic 
mean-field density |^. Figure 11a) shows the true density distribution of "^^Ni together with the 
distribution fitted by the use of hfteen {M = 15) Kamimura-Gauss basis functions. They agree each 
other quite closely, as they could not be identified separately by eye. = 4.181 [fm] is obtained 
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from these density distributions. 

We use the same size parameters, ri = 0.7 [fm], and tm = 5.6 [fm], as those that we used for ®®Ni 
with M = 15 (Table 1). The completeness errors for M = 10 and 15 are shown in Fig. 11 b), which 
are very close to those for ®®Ni. Therefore, once we determine the size parameters for a given M for 
a known density distribution, such as that for stable nuclei, we could safely use them for the isotopes 
of the nucleus. 

The numerical results of the least square fittings are shown in Fig. 12 a) and b). We use i? = 10 [fm] 
and M = 10. The size parameters, ri = 1.0 [fm], and tm = 5.5 [fm], are the same values as those of 
®®Ni (Table 1). The pseudo data are fitted well, especially in the forward direction, as seen in Fig. 12 
a). Figure 12 b) shows that the fitted density distributions fluctuate in the central region, but are 
stable in the surface region. Unfortunately, obtained from the fittings is slightly smaller than 

that of the true distribution shown in Fig. 11a), but within the fluctuation. 

Similarly to the discussions in the previous subsection, we also show the statistical uncertainty 
and the systematic error in Figs. 13 and 14. Figure 13 confirms that the density distribution for 
this neutron-rich nucleus can be determined well in the surface region, which is consistent with the 
result obtained by using the Fourier-Bessel basis functions (Fig. 15). M = 7 is the best choice for this 
pseudo-data-set group for the Fourier-Bessel basis functions of i? = 10 [fm]. 

One should be careful for the statistical uncertainty in the central region. By comparing Fig. 13 
with Fig. 15, one can see that the statistical uncertainty of the Kamimura-Gauss basis functions in 
the central region fluctuates much larger, about an order of magnitude, than the Fourier-Bessel ones. 
As we mentioned in Introduction, since the statistical uncertainty should carry the information only 
on the quality of data, any basis function has to show the statistical uncertainty distribution of similar 
quality for the same data. The error estimation of the Kamimura-Gauss basis functions in the central 
region should be treated with caution. 

The behavior of these uncertainty and error is similar to that for the stable nucleus examined in 
the previous subsection. We thus expect that the surface region of nuclei could be better determined 
than the central region for relatively low intensity beam. 

Before closing this section, we should add a comment. To discuss how the density distributions or 
the radii change in the isotope chain, the difference between the one from the neighboring nucleus is a 
good measure, because the model dependence could be canceled by the subtraction procedure js^, 

In this work we are interested in how precisely we can determine the density distribution itself of the 
neutron-rich nuclei with low intensity beam. We leave the problem concerning the discussion on such 
differences for our future work. 
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5 Summary and Discussion 

The matter density distribution in the surface region is shown to be determined well by a relatively 
low intensity beam in the order of 10 [sec“^], which is an expected intensity of the ^®Ni beam at the 
RIBF of RIKEN. 

We have demonstrated how the determination of the density distributions becomes poorer, espe¬ 
cially in the central region, as the quality of data becomes worse. The demonstration is made through 
a simulation by the use of three pseudo-data-set groups of proton-nucleus elastic scattering, corre¬ 
sponding to beam intensities expected to be available at the RIBF. Our discussion here is qualitatively 
in accord with what has been found in electron scatterings 1^ 1^, |ll|, |l^, 1^, |l^ , and also with what 
has been done in proton scatterings ©il- 

The matter density distribution of neutron-rich nuclei is determined by fixing the parameters in 
the basis functions to those of the neighboring stable nucleus. This procedure is expected to be a 
practical way of treating neutron-rich nuclei, as the beam intensities of neutron-rich nuclei tend to be 
less and the determination of their matter density distributions will be poorer. We have demonstrated 
this procedure to be feasible by taking ^®Ni as an example. 

In this work, we have focused on the question, how well the nuclear density distribution is deter¬ 
mined from experimental data with the quality expected at the upcoming facilities, assuming that 
the proton-nucleus potential is known with no uncertainty. The uncertainty for this must be included 
in the practical evaluation of the uncertainty in the determination of the matter density distribution. 
The question of the potential will require a close examination beyond the use of t-p form, especially 
as the incident energy becomes lower. Also, a special attention may have to be given to neutron-rich 
nuclei 

We expect that the uncertainty associated with the use of the eikonal form instead of the partial 
wave decomposition is small, though the absolute magnitude itself could be substantial. In order to 
clarify this and the uncertainty associated with the potential, they will be the issues that we wish to 
address in future work. 
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APPENDIX 


A Kamimura-Gauss Basis functions 


Here we summarize the formulae associated with the Kamimura-Gauss basis functions, which have 
been used by Kamimura in various applications, for example, in few-body calculations 
The basis functions are expressed as 


fn{r) =exp{-r'^/rl) 9{R-r), (n = 


(A.l) 


where the size parameter, r„, is determined through a geometrical progression, r„ = ria^~^ (n = 
1, • • •, M). 6{x) is the step function, and i? is a scale parameter constrained by the kinematical con¬ 
dition |l^. Note that the standard version of the Kamimura-Gauss basis functions has no restriction 

of tm < i? [|3- 

The inner product of two basis functions is dehned by 


where 


nOQ 

{fn,fn+m)i = 47r/ r'dr/„(r)/„+m(r) 

Jo 


2tt 


d+i)/2 

G^n.m 


{T ((/ + l)/2) -T{il + l)/2, 


^n,m 


1 1 
^ l2 

n n-\-m 



(A.2) 

(A.3) 


(A.4) 


and I is the power of the phase volume and is set to be a non-negative integer. r( 2 :) (= r(z, 0)) is the 
Gamma function, and T{z,p) is the incomplete Gamma function of the second kind. 

When the normalization is suitably chosen, the expression of the inner product becomes simple. 
If ifn, fn)i=2 = 1 for fnir) = fn{r), we obtain 


= (r(3/2) - r(3/2, 2i?Vr2)) . 

Note that when R oo, the expressions become greatly simplified: 

1 _^3/2 ^ 

~ 2 V 2 


and 


fn-\-m)l—2 


2V2a3™/2 
(1-f a2'")3/2 ■ 


(A.5) 


(A.6) 


(A.7) 


The basis functions have the following properties: 1) They are not orthogonal, but are linearly 
independent. One can easily verify the linear independence by demonstrating the Wronskian, to be 
non-zero except at the infinity. 2) The derivative with respect to r vanishes at r = 0, and thus the 
derivative of the density expanded in terms of the basis functions is also zero at r = 0. Note that this 
property is the same as that of the Gaussian basis functions used by Sick |Q. 3) The coefficients of the 
Kamimura-Gauss basis functions can be negative, while Sick restricted the coefficients of his Gaussian 
basis functions to be positive . Actually, the coefficients in our case tend to appear almost equally 
in both signs. 4) The set of the Kamimura-Gauss basis functions does not form a complete set, even 
if the inhnite number of terms are used. Note that a finite number of terms is always used in practice, 
and the issue of completeness is an academic question. 5) They depend on two parameters, such as 
Ti and tm- 
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In our application, we choose ri and tm by minimizing D{Cn) of Eq. (A^) below. Though ri 
depends strongly on M, it is found to vary little for different isotopes. The expansion coefficients for 
the normalized basis functions, C„, are determined by minimizing 


D{Cn) = IPtrueW -PmWI 


M 


= \PtTne{r) -'^Cn fn{r) 


n—1 


that is, by satisfying dD{Cm)/dCm = 0. This condition leads to 


M 


fmir) Ptrue(r) = X! {n,m = 


We thus obtain 


\n—l 


M 




(A.8) 

(A.9) 


(A.IO) 


(A.ll) 


where 


, 1 ,™ = 47r / r'^dr /„(r)/m(r). 

Jo 

PCG 

Vn = in r’-dr f,n{r) PtTue{r). 

Jo 


(A.12) 

(A.13) 


Note that one can also follow the same procedure in the case of the coefficients, C„, for non-normalized 
basis functions. The phase-volume parameter, I, is chosen by depending on how much weight one 
wishes to impose on the out-region. In the present analysis, we choose Z = 0 in order to include the 
information of the whole region equally. 

Finally, Eq. (A.ll) is inverted numerically to yield 


M 


Cti — Pm,n 


(A.14) 


n—1 
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B Eikonal Expressions 

For the purpose of completeness and clarification, we summarize the eikonal expressions [IT^ , which 
we have used. 

The differential cross section of the proton-nucleus elastic scattering in the center-of-mass frame 
is expressed as 

^ ITMI^ (B.l) 

with the kinematical factor Flcin = |pP/(271')^. q (= p — pO is the three-momentum transfer, where p 
and p' are the momenta of the incident and of the outgoing proton, respectively. r(q) is the t-matrix 
expressed in the eikonal form, 


T(q) = y db exp{tq • b} (1 - exp{ix(b)}), (B.2) 

nOC 

= 2tt bdb Jo{qb) (1 - exp{ix(&)}), (B.3) 

Jo 

where Jo{qb) is the zero-th order Bessel function. The profile function, x(b), is given by 

27r 

X(b) = T-T/(O) / dz p{h,z'), (B.4) 

where p(r) is the one-body matter density distribution, and r = (b, z). b is the impact parameter 
vector, b _L p, and the integration in z is taken in the direction of the incident momentum, p. /(O) 
is the nucleon-nucleon scattering amplitude in the forward direction, q = 0. 


/(O) = ^(* + a) 4 r- 


(B.5) 


The expression of x(b) above corresponds to the first-order optical potential in the Glauber approx¬ 
imation. The Coulomb interaction is not included here, but could be included straightforwardly 

EH- 


For the numerical calculation of Fig. 1, we use the parameters of Eq. (B_^) as follows [p2[: 


IpzabI = \/{K + m)2 - rn? = 1.75 [GeV/c]. 

a = -0.0687, = 43.1 [mb]. 


(B.6) 

(B.7) 


where m is the nucleon mass, and K is the kinetic energy of the incident proton, which is K = 
1.047 [GeVj. We use these parameters for all the numerical calculations in this work. 

We calculate the cross section of proton-®®Ni elastic scattering as an example in Subsection 3.1 
and Section |^. We take the three-parameter Fermi distribution, Eq. (O), for the matter density 
distribution of ®®Ni to create the pseudo data in Subsection 4.2, because the experimental data of 


differential cross sections of electron-nucleus elastic scattering for various nuclei are fitted well with 
this form, including ®®Ni |Q, and because the matter density distribution of ®®Ni is considered to be 
quite similar to the charge density distribution. Actually, the experimental data of differential cross 
section of proton-®®Ni elastic scattering are fitted well with this form ]3l[ |. 

One should be careful for large-r behavior of the three-parameter Fermi distribution, Eq. 
because it becomes negative for r > Tc if w is negative. That is, p{rc) = 0 1 + w{rc/c)'^ = 0, for 

= cj^/—w. For our choice of the parameters, Cc = 11.91 5 [fm ]. We put p{r) = 0 for r > Tc- 

With the fitted density distribution, p%{r), of Eq. (^) in terms of the Kamimura-Gauss basi s 
functions defined by Eq. ( t3.4|) , we can perform the integration of the profile function in Eq. ( |B.4| ) 
analytically due to the nature of the Gaussian functions: 


27r 


M 


^ exp{-bVc^} 


n=l 


X {r(l/2)-r(l/2,(z,nax/r„f)}, 


(B.8) 
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where r(l/2) = and — b^. T{z,p) is the incomplete Gamma function of the second 

kind. This integration reduces our numerical task for the cross sections, and brings us the numerically 
stable results. 


18 



R = 

10 [fm] 

M 

Am 

5 

58.115 

10 

57.974 

15 

58.001 

20 

57.999 


Table 3: Dependence of the integrated value of the density Am of the Fourier-Bessel basis functions 
on the number of the basis functions M. The radial range of the fit is i? = 10 [fm]. 


C Fourier-Bessel Basis Functions 


Here we summarize some formulae concerning the Fourier-Bessel basis functions, and show its com¬ 
pleteness error introduced in Subsection 3.3.1. The Fourier-Bessel basis functions are defined by 0 


fn{r) = joimrr/R) 0(i? - r), 


(C.l) 


where jo{x) is the 0-th order spherical Bessel function, which is expressed by joix) = sinx/x. i? is a 
scale parameter determined by the discussion under eq. (O). 0(x) is the step function. 

If we define the inner product as 


(/n,/m);=2 = 47r / r'^dr fn{r)fjn{r), 

Jo 

we obtain the normalized version of this basis function as, 

/n(r) = Nn fn{r), 


(C.2) 


(C.3) 


where Nn = ^n^'K /{2R^). 

Since the Fourier-Bessel basis functions are orthogonal to each other, i.e., (/„, fm) = Sn,m- we can 
easily obtain Cn by the following equation: 


Cn = 47r 


r'^dr fn{r) ptrue{r). 


(C.4) 


As in Subsection O, we show the numerical results of Ap^™(r), defined by Eq. ( 3.10D , for the 
Fourier-Bessel basis functions of i? = 10 [fm] with the number of terms, M, being 5, 10, 15, and 20, in 
Fig. 16 a). This figure corresponds to Fig. 2 a). Here we also systematically determine the expansion 
coefficients of the Fourier-Bessel basis functions, Cn, of Eq. (C.4). Roughly speaking, the magnitude 
of the coefficients becomes smaller as n increases, because the overlap of /n(r) with ptrue(x) becomes 
smaller due to the oscillatory nature of the basis functions as one can see from Eq. (C.4). 

Although each term of the Fourier-Bessel basis functions and the Kamimura-Gauss ones shows 
different behavior in r, it is interesting to see that Ap^™(r) of the both bases show similar oscillating 
patterns in r. The convergence at each point becomes better as M increases for the both cases, and 
the convergence for the Fourier-Bessel basis functions become better beyond M = 15 in contrast to the 
Kamimura-Gauss ones. Up to M = 15, the both bases have almost the same quality of convergence. 
The uncertainties for M > 10 are the order of 10“^ or less, and for M > 15 are the order of 10“^ or 
less in the region of 1 [fm] < r < 5 [fm], which covers the surface region of nuclei. The uncertainty 
becomes worse in the region of r > 5 [fm], but it is still the order of 10“^. 

As in the completeness calculation in Subsection 4.1, we impose no normalization condition (t5.2|) 


on pM^r), Eq. (3.11), here, and thus the integrated value of pm^x) also serves as a measure of goodness 


of PM{r). As in Table 1 in Subsection 4.1, Table 3 lists the integrated value. Am- Note that the exact 


integrated value is 58 in this case, too. As one can see from this table, Am stays in about 0.2 %- 
fluctuation for all M. They are all tolerable. For M larger than 5, Am can be determined within the 
error of 0.05 %. 
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Finally, we show the numerical results of the corresponding cross sections for each M in Fig. 16 b). 
This figure corresponds to Fig. 3. The Fourier-Bessel basis functions give satisfactory convergence for 
M > 10, similar to the Kamimura-Gauss ones, and, if we restrict ourselves to see 9 < 10°, we already 
obtain good convergence even for M = 5. The solid curve (M = 20) is almost on the dotted curve 
(drawn by the original ptrue(^)) for < 20°, but unfortunately the difi'raction pattern is broken for 
0cm > 20°. The other curves, M < 20, are also differs from the dotted curve in this large-angle region. 

Such broken patterns of the Fourier-Bessel basis functions could be numerical artifacts coming from 
the numerical integration, which is almost inevitable for numerical fittings, while they do not appear 
for the Kamimura-Gauss basis functions. Thus, when we apply the Fourier-Bessel basis functions to 
those fittings, we should be careful for cross sections of the large scattering angle. They most probably 
affect the central region of the density distributions. For 9 < 20° and M > 10, both the Fourier-Bessel 
and Kamimura-Gauss basis functions give good convergence. 
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D Small Scattering Angles 


Here we explain how to reduce the separation of pseudo data in Oq, from 1.0° near the forward direction. 
The prescription to create pseudo data is written in Subsection 3.2. To make the way to shrink the 
step of pseudo data systematic, we define the minimum uncertainty, Cmin- Since we are to assign 
a statistical uncertainty more than this, we restrict < l/fVpD(0a) at each data poi nt, where 
iVpD(0ct) is the yield count of the pseudo data at 0 q, defined in the same way as in Eq. (^). We 
obtain A0 q by keeping 


A0O “2- ^ A0q N'pji{9a)- (D.l) 

^min 

In the region where the theoretical cross section changes drastically this formula overestimates or 
underestimates the step. 

For large A^pD(^a), A0„ becomes smaller than a bin of the detector. In order to avoid such an 
unrealistic estimate, we define the minimum step, A^mirn and make A0 q not less than that. The 
minimum step is expected to be 0.3 [mrad] (= 17.2° xl0“^) in the proton-fixed frame at the RIBF. 
Since our cross sections are calculated in the center-of-mass frame, we are careful for the inverse 
kinematics of the RIBF. We take A0i„in = 2 x 0.3 [mrad], which is safe if we transform A^min in the 
proton-hxed frame to the center-of-mass frame. The requirement to identify a proton-nucleus elastic 
scattering makes the minimum step very small. In the forward direction 6 can be determined precisely 
at the RIBF. Although the data have some uncertainty in 9 in the large 0cm region, we neglect it in 
this work. 
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E Distribution of the Fitted Results 


According to the statistics, the least square, of Eq. ( ll.3[) is known to obey the x^-distribution 
defined by 


J 2 ^ ^exp(-xV2), (x^ > 0), 

I 0, (x^ < 0), 


where v (> 0) is the degree of freedom, provided that each term consisting the least square of Eq. (3.3), 
(^p.data/^^ j obeys the standard normal distribution Q. 

We use this fact to examine the consistency of our simulation by calculating the following quantity. 



dv f^{v). 


(E.l) 


We compare it with our numerical results in Table 4. 

We take the pseudo-data-set group B as an example. The number of data points of each pseudo¬ 
data set is 40, 41, and 42. Such a fluctuation of the number of data point occurs, because some pseudo 
data are rejected from the pseudo-data set when they become less than 10“^ [mb/str]. We use the 
pseudo-data sets of 41 data points, because this is the largest portion of the group. The degree of 
freedom becomes 31, as M = 10. The case contains 15 pseudo-data sets out of 25. 

Table 4 shows that our fitted results distribute consistently with the calculated x^-distribution, 
and that our simulations are in accord with statistical considerations. 


Xa 

Fixt) calc. 

F{xt) PDSG-B 

20 

0.064 

0.067 (= 1/15) 

30 

0.48 

0.33 (= 5/15) 

40 

0.87 

0.87 (= 13/15) 

50 

0.98 

1.00 (= 15/15) 


Table 4: Comparison of calculated values of integrated x^ distribution with the corresponding quantity 
of the pseudo-data-set group B. The degree of freedom is 31. The bracket in the right most column 
shows the fraction which satisfies the condition. 
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Figure 1: Comparison of the theoretical cross section with the real experimental data. The solid curve 
is the theoretical cross section of the proton-®®Ni elastic scattering without the Coulomb interaction. 
The crosses with bar are the experimental data of SATURN at Saclay. The kinetic energy of the 
incident proton is 1.047 [GeV]. 
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Fig. 2a) 
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Fig. 2b) 



r [fm] 


Figure 2: M-dependence of the completeness error of the Kamimura-Gauss (KG) basis functions 
a) with the best-fit parameters, b) with a crude choice of ri and tm, for M = 5, 10, 15, and 20. R 
= 10 [fm]. The dotted curve is for M = 5, the dashed one is for M = 10, the dash-dotted one is for 
M = 15, and the solid one is for M = 20. 
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Fig. 3 



Figure 3: M-dependence of the cross section with the Kamimura-Gauss (KG) basis functions for 
M = 5, 10, 15, and 20. i? = 10 [fm]. The dotted curve is for M = 5, the dashed one is for M = 10, 
the dash-dotted one is for M = 15, and the solid one is for M = 20. Another dotted curve, as it is 
quite difficult to be seen, is obtained by using the original ptrue(?')- 
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Figure 4: An example of pseudo-data set of the pseudo-data-set group (PDSG) A. Those simulate 
a series of measurements of the cross section of proton-^®Ni elastic scattering. The solid curve is 
obtained by using the original ptmeir). 
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Fig. 5a) 



Fig. 5b) 
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Figure 5: Results of the least-square fitting for the pseudo-data-set group (PDSG) A of ®®Ni. R 
= 10 [fm] and M = 10. a) 25 fitted cross sections (solid curve) and all the pseudo-data (cross with 
bar), b) 25 fitted density distributions (solid curve), and the original distribution (dotted curve). The 
density distributions are drawn in the form of p(c). 
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Fig. 6b) 


b) Fitted Density Distributions 
(R = 10 [fm]: M = 10) 


[tm]: 3.781 ± 0.018 


Figure 6: Results of the least-square fitting for the pseudo-data-set group (PDSG) B of ®®Ni. R 
= 10 [fm] and M = 10. a) and b) are the same as Fig. 5. 
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Fig. 7b) 



Figure 7: Results of the least-square fitting for the pseudo-data-set group (PDSG) C of ®®Ni. R 
= 10 [fm] and M = 10. a) and b) are the same as Fig. 5. 






Figure 8: Results of the least-square fitting for the pseudo-data-set group (PDSG) C’ of ®®Ni. R 


= 10 [fm] and M = 10.. a) and b) are the same as Fig. 5. 






Fig. 9a) 
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Figure 9: a) Statistical error, and b) systematic error of the fitted density distributions for ®®Ni. We 
divide them by ptrue{r) at each r. The solid curve is for the pseudo-data-set group (PDSG) A, the 
dotted one is for the group B, the dashed one is for the group C, the dash-dotted one is for the group 
C’. 
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•e 10: a) Statistical error, and b) systematic error of the fitted 
IS are the same as in Fig. 9. 
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Figure 11: a) Matter density distribution of ^®Ni. The dotted curve is the original density distribution 
obtained in the relativistic mean-field theory, i? = 10 [fm] and M = 15. The dashes curve is obtained 
by fitting of the Kamimura-Gauss (KG) basis functions, b) Completeness error of ^®Ni for M = 10 
and 15. The solid curve is for M = 10, ri = 1.0 [fm] and tm = 5.5 [fm] (Table 1). The dashed curve 
is the case of M = 15, n = 0.7 [fm] and tm = 5.6 [fm] (Table 1). 
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Figure 16: M-dependence of a) the completeness error of the Fourier-Bessel (FB) basis functions and 
of b) the cross section with the same basis functions for M = 5, 10, 15, and 20. i? = 10 [fm]. The 
dotted curve is for M = 5, the dashed one is for M = 10, the dash-dotted one is for M = 15, and the 
solid one is for M = 20. Another dotted curve in b), as it is quite difficult to be seen, is obtained by 
using the original /9true(r). The curves are the same as in Figs. 2 and 3. 
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